Velocity and Velocity-Difference Distributions in Burgers Turbulence 
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We consider the one-dimensional Burgers equation randomly stirred at large scales by a Gaus- 
sian short-time correlated force. Using the method of dissipative anomalies, we obtain velocity 
and velocity-difference probability density functions, and confirm the results with high-resolution 
numerical simulations. 
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1. Introduction. The one-dimensional Burgers equa- 
tion with a random external force, 



tit + uu x = vu xx + 



f(x,t), 



(1) 



is a simple model of turbulence, where u(x, t) is a one- 
dimensional velocity field, and v is small viscosity. The 
external force is assumed to be Gaussian, with zero mean 
and white-in-time covariance, 



(f(x,t)f(x',t')) = K(x-x')S(t-t'), 



(2) 



where n(y) should be specified. We will assume that it 
is an analytic function at y — > with the characteristic 
scale L; non-analytic forms of k have been considered as 
well Q . This model has attracted considerable attention 
since it is believed to be exactly solvable El II I! 0, IE IE 

II II IE EE El El E! El. 

Let us introduce the characteristic function of the Ap- 
point velocity distribution, Zn(Xi,xi; \N,XN',t) = 
(exp(Aiu(a;i, £) + •.. + Xnu(xn, £))}. As was shown in Q, 
it obeys the master equation that one derives by differ- 
entiating this function with respect to t, and by using 
Eq. Q, 



dZ 2 
dt 



d 1 dZ 2 1 



Xj )Z 2 + D, (3) 



where the summation over repeated indices is assumed. 

The last term in Eq. JSJ denotes the contribution of 
the dissipative terms. Although the method of our pa- 
per is valid for the general case of iV-point character- 
istic function, we will concentrate on the case N = 2, 
where the solution is easy to find. In this case, D = 
v{(Xiv, XlXl + \2U X2X2 )exp(\ 1 u(x 1 ,t) + \ 2 u(x 2 ,t))). This 
term does not vanish in the limit of infinitely large 
Reynolds number, v — ► 0; rather, it has a finite value, 
because the velocity field develops singularities, shocks, 
where large velocity gradients are balanced by the small 
viscosity. It has been proposed in Q that the origin of 
this term is analogous to the origin of anomalies in quan- 
tum field theories, and it has been suggested that this 
term can be expressed linearly through the Z function. 
Namely, the following ansatz should be true inside any 



N > 2 point correlation function 



A = lim v\u xx e Xu{x) = 

I/->0 



a 6-1 d 
2 



A dx 



,A«(x) 



(4) 



The right hand side of |Q| contains the simplest linear 
in e Xu terms that are consistent with translation, scale, 
and Galilean invariance H IE]. The higher-derivative 
terms are not allowed since they would change the struc- 
ture of Eq. Q and may lead to additional, non-physical, 
solutions. The parameters a and b, which we will call 
'anomalies,' should be found from the requirement that 
the stationary solution of Eq. © correspond to a posi- 
tive, finite, and normalizable PDF, similar to the eigen- 
value problem in quantum mechanics. The other condi- 
tion is non-positivity of the dissipation, D < 0. These 
two conditions restrict the values of the anomalies con- 
siderably. It has been shown in that the admissible so- 
lution exists for a one-parameter family {6, a(&)}, where 
3/4 < b < 1 and the corresponding interval for the a 
anomaly is > a > —0.45. 

The corresponding scale-invariant solution for the ve- 
locity difference PDF, w(Au/y), where Am = u(x\) — 



u(x 2 ) < u r 



and 



y 



X\ 



i2 « i can be found 



from ©• It has a peculiar structure that was first 
qualitatively established in numerical simulations 
and was confirmed by a variety of analytical meth- 
ods H H 0, IE IE ■ It decays hyper-exponentially fast 
for the large positive argument, w(z) oc exp(— z 3 /3), and 
has an algebraic tail for the large negative argument, 
w(z) oc |z|~ 2fc_1 [see section 2]. This behavior is under- 
stood since due to the nonlinear term in Eq. 0, positive 
gradients decrease in the course of time, while negative 
gradients become stronger. Therefore, large negative ve- 
locity gradients are more probable than large positive 
ones. 

According to the allowed values of the anomaly b in 
model |4jl. the distribution density should decay slower 
than |^| 3 for z — > — oo. However, the geometric consid- 
eration of pre-shock events producing large negative ve- 
locity gradients, proposed in [21 1^1, suggested the asymp- 
totic |z|~ 7 / 2 . Lagrangian simulations of o nly the left tail 
of the PDF agreed with the latter result [HE! El- 
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So far, no analytical derivation of the full velocity- 
difference PDF decaying as |z|~ 7 / 2 was available. The 
direct numerical simulations of Eq. QJ, carried out 
in 0, , have not convincingly reproduced such a PDF 
either. The form of the PDF has therefore remained a 
subject of controversy. 

In this Letter we propose that the PDF is not unique. 
In an infinite system, where the Galilean invariancc 
holds, the anomaly has the form (J2J. In the finite-size sys- 
tems, investigated in 0,EJ0|) the global Galilean in- 
variancc is broken and the weak Galilean principle should 
be applied. In this case the form of the anomaly is dif- 
ferent from 10}, and the solution with the asymptotics 
z|~ 7 / 2 may be allowed. We derive the corresponding ve- 
locity and velocity-difference PDFs. Next, we conduct 
extensive numerical simulations of Eq. , using a high- 
resolution shock capturing scheme. The results agree well 
with our analytical prediction. Section 2 presents the 
theory, section 3 - numerical simulations. 

2. Velocity and velocity-difference PDF's. The resolu- 
tion of the controversy mentioned above is based on the 
fact that the form Q of the anomaly is dictated by the 
strong form of the Galilean principle (see below) while the 
arguments of as well as numerical results |12lI13lI14I | 
are correct only when the weak G-principle applies. The 
strong Galilean principle is the requirement that in the 
infinite system, taking velocity to zero at infinity does 
not break the Galilean symmetry in the middle. This 
requirement is analogous to the absence of the sponta- 
neous symmetry breaking in ma gnetic sy stems. For the 
systems considered in [J 0, H2I Il3l Hif this strong G- 
principle obviously does not hold. The reason is that 
the size of the systems is of the same order as the cor- 
relation length of the force. Therefore, the velocities in 
the middle of the system are correlated to those at the 
boundary. The same argument works for the periodic 
case, in which the Galilean symmetry is broken by the 
condition Uq = J udx — 0. 

So, the results in @ and 0IEE1I3 refer to 
the different dynamical systems and that explains the 
discrepancy. The natural question is whether we can ap- 
ply the methods of [2j for the case in which the strong 
G-symmetry is broken. This is important because the 
method of dissipative anomalies is the only known non- 
perturbative way to solve Burgers turbulence. In this 
case we still have the weak G-symmetry, which says that 
all correlation functions are G-symmetric, provided that 
we average them over uq. This is of course always true, 
but how to impose this condition? 

Let us notice that the correlation functions calculated 
at the fixed uq have the property 



,t)\ 



(5) 



the form: 

(e£ = 2n5 (J2 X i) (e^'^'^o- (6) 

The master equation must allow this ansatz, namely 
the proportionality to the 5 (Y) X-j) must be consis- 
tent with the equation. This is obviously so with the 
anomaly (@J, but this also allows a more general form of 
the anomaly. Let us look at the expression 



A = A + cX-^-e > 

dX 



(7) 



it adds the term D = c^XjdZ / dXj to the right hand 
side of Eq. ©. When Z = 6 A,-) Z , we have D = 
c (J2 XjdZo/dXj — Zq) S (J2 Aj), where we used the iden- 
tity xS'(x) = —S(x). It is important to notice that any 
other function of A or higher power of d/dX would break 
the ansatz. Anomalies of the form (JJJ) have been con- 
sidered in relation to incompressible turbulence by V. 
Yakhot [private communication]. 

The restrictions on the coefficients a, 6, and c can 
be obtained from the equation for the velocity-difference 
PDF. To write this equation, substitute the dissipative 
term Q into Eq. ©. In the Galilean-invariant limit, the 
statistics of velocity differences are separated from those 
of the mean velocity field. Let us change the variables, 
Ai = A + ^i, A 2 = A — /i, xi = X + y/2, x 2 = X — y/2, and 
assume that A <C /J, and y <C L. The latter condition 
allows us to expand n(y) — k,q — K 2 y 2 /2. In this limit, 
the master equation becomes 



d 2 Z 2 2b dZ 2 
d[idy n dy 



- 2A 2 n Z 2 - 



2 2 

y 



Cfl- 



9£2 



cA 



dZ 2 
~dA 



aZ 2 



(8) 



where we consider the stationary case, dtZ 2 = 0. In the 
finite-size system with uq — 0, we assume the weak G- 
symmetry and look for the solution in the form Z 2 = 
Z+(A)Z-(/j,,y). [The strong G-symmetry would imply 
Z+(A) ~ (5(A) and the solution might be different.] Sub- 
stituting this into Eq. (JSJ we get the following set of equa- 
tions for the characteristic functions and Z_ : 



2A 2 n Z + = cA 



dZ 4 



(9) 



2b dZ- 



n 2l x 2 y 2 



d\xdy [i dy 



dA ' 

dZ- , s 

c/j— — + aZ_. (10 



The solution of Eq. (|5J corresponds to a positive and 
normalizable PDF only when c < 0, and the solution is a 
Gaussian. Its Laplace transform then gives the velocity 
probability density function 



Hence, if we integrate over uq the correlations must have 



P {l ; ]= /^iexpj^- 

7TKQ V K 



(11) 
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where U = [u(xi) + u(x2)]/2. Since the point separation 
is much smaller than the force correlation length, P(U) 
becomes the one-point probability density function of the 
velocity field. Strictly speaking, expression ljTT|) is phe- 
nomenological, it is expected to match the velocity PDF 
only in the region U <C U rms ; for larger U the velocity 
PDF may be not universal. We compare this result with 
the simulations in section 3. 

Equation QlOj l is more complicated, but it can be sim- 
plified if we are looking for the solution in the scale- 
invariant form, Z-(p,,y) = where x = fiy, 

x<&" + (1 - 2bW - — x 2 <$> = a<f> + cx— $. (12) 
2 ox 

To rewrite this equation in the velocity space, let us 
Laplace transform the $ function, 

*(x)=J W (z)eM*Z)dz, (13) 

— oo 

The function w(z) is then related to the velocity- 
difference probability density function, W(Au, y), as 
W(Au,y) — w(Au/y)/y. Equation for w is readily ob- 
tained by substituting l|13(l into |1 2jl. Introducing the 
dimensionless variable z = z(k,2/2)~ 1 ' s , we thus get 

w" + z 2 w' + (1 + 2b)zw = -aw + c(zw)' , (14) 

where the derivatives are with respect to z, and a = 
a(K 2 /2) -1 / 3 and c = c(k 2 /2)~ 1 / 3 . Asymptotics of the 
solution at \z\ — > oo can be found from Eq. I|14f> . 

w ~ z (2fc " 1) exp(-z 3 /3 + cz 2 /2), z -> +oo, (15) 
w ~ |z|-( 2b+1 \ z^-oo, (16) 

they have the form that we already discussed in the in- 
troduction. 

Ideally, we should be able to derive all the anomalies 
a, b, and c from our theory. At present we are not able 
to do this, so in the next section we turn to numerical 
simulations. 

In general, the anomalies a, b and c should satisfy the 
requirement that Eq. (|14|) have positive and normalizablc 
solution, and that the dissipation is non-positive, aw + 
2(b-l)zw-c(zw)' < 0. This eigenvalue problem thus 
has two-parameter family of solutions. 

Note that the geometric approach of 0, H, El 
Q] , proposed the asymptotic for negative velocity differ- 
ences, w ~ |z| -7 / 2 . This condition is satisfied by b = 5/4. 
This value was not allowed by the condition of non- 
positivity of the dissipation in the theory based on 10}, 
but with the new expression for the anomaly Q this 
value is admissible. For b — 5/4, Eq, l|14|l can be solved 
numerically, leading to the approximate window of ac- 
ceptable values, —0.37 < c < —0.025, which turns out to 



be consistent with (|llfl . We now find the values of the 
anomalies from numerical simulations of Eq. Q . 

3. Numerical results. We perform numerical simula- 
tions of Eq. P in a periodic interval of length L = 1 
containing 10 6 grid points. The used numerical scheme 
is the standard shock capturing, total variation dimin- 
ishing (TVD) scheme with non-linear, limited solution 
reconstruction [l5j |. For non-linear, scalar equations the 
scheme is provably TV-stable, monotonicity preserving 
and second-order accurate in the L p , p < oo, norm. To 
resolve discontinuous solutions the scheme uses nonlinear 
numerical viscosity that appears only inside the shocks. 
This dissipation mechanism ensures that the shocks are 
restricted to 2-3 grid points. 

The external force is generated as 

m 

f(x) = V] A(k) [£ fc sin(27rfcr) + % cos(27rfcr)] , (17) 

k=l 

where and rjk are independent Gaussian random vari- 
ables with mean and covariance 1. The amplitude of 

1/2 

the force is A(k) — k [exp(— /c 2 /n 2 )/r] , where r is the 
force renewal time; in our simulations we had t = 10~ 5 , 
while the integration time step was to = 10~ 7 . We used 
m = 10, and performed a series of runs for n — 5, 10. In 
this paper we mainly present the results averaged over 5 
independent runs with n = 5, and with the integration 
time t = 2 each. 

We found that the one-point PDF can be best matched 
by the solution of Eq. ((lip, with the choice c = —0.3. 
Formula (|TT1) gives 

„ /rrN /o.3 4 /3 / 0.3K2 /3 Tr2 \ 

and Fig. shows its good agreement with the simula- 
tions for both n = 5 and n = 10. As expected, the best 
agreement is seen in the top parts of the curves; the tails 
of the PDF seem to be non-universal. 

Once the c anomaly is given, we find the anomaly b 
by matching the numerically obtained velocity-difference 
PDF with the solution of Eq. JUj), see Fig. 10 El The 
best match to the whole PDF is given by b pa 1.19 
which is slightly less than the value b = 5/4, predicted 
in and leads to the left-tail exponent 

-3.38. The solution of Eq. JT3J with the value b = 5/4 is 
less consistent with our simulations. Although the reason 
for this discrepancy is quite intriguing and is not clear to 
us, our numerical resolution is not high enough to dis- 
tinguish between the two asymptotics, —3.5 and —3.38. 
Therefore, more extensive numerical work is required to 
understand this difference. 

Once the anomalies c and b are found, the a anomaly 
is obtained from the eigenvalue problem (|14f) . which 
gives a w —0.65. 
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ONE-POINT PDFS 



VELOCITY-DIFFERENCE AND VELOCITY-GRADIENT PDFS 





FIG. 1: Numerical and analytical one-point PDFs for two 
runs, n — 5 and n = 10. The broader curves correspond 
to the run n = 10. The plot is in log-lin scale. 



VELOCITY-GRADIENT PDF 




FIG. 2: Velocity-gradient PDF for the run n = 5. The 
analytical curves correspond to the solution of Eq. 11411 
with b = 1.19 (dashed line) and b = 5/4 (dotted line). In 
both cases, c = —0.3. The plot is in log-lin scale. 



4 Conclusions. We have found that our analytical 
model (|14fl with parameters a « —0.65, b w 1.19, and 
c ~ —0.3 provides an excellent fit to the numerically 
obtained velocity and velocity-difference PDF's of Burg- 
ers turbulence. As the important numerical check, one 
can directly find the form of the anomalies by numeri- 
cally constructing conditional probabilities of shock am- 
plitudes for given velocity gradients on the left and on the 
right of the shock Q . We plan to address this question 
in the future with more extensive simulations. 

We also believe that the method of dissipative anoma- 
lies can be generalized for turbulence with pressure and 
can be applicable in more than 1 dimension. We cur- 
rently conduct the corresponding simulations, the results 
will be reported elswere. Practically, the Burgers model 
is expected to work for strongly compressible flows, i.e., 
in astrophysical application, such as supersonic turbu- 



FIG. 3: The left tails of the velocity-gradient PDF (lower 
curve) and of the velocity-difference PDF (upper curve, 
point separation Ax = 10~ 3 ). For reader's orientation, 
the straight line has the slope —3.4. 



lence in cold molecular clouds [lil Il7| . 
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